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Abstract 

In this paper we consider the matching coefficient of the vector current between 
Quantum Chromodynamics (QCD) and Non-Relativistic QCD (NRQCD) to three- 
loop order in perturbation theory. We evaluate the fermionic corrections containing 
a closed massless fermion loop. The results are building blocks both for the bottom 
and top quark system at threshold. We explain in detail the methods used for the 
evaluation of the Feynman diagrams, classify the occurring master integrals and 
provide results for the latter. The numerical effects are significant. They have the 
tendency to improve the behaviour of the perturbative series — both for the bottom 
and top quark system. 

PACS numbers: 12.38.Bx 14.65.Ha 



1 Introduction 



One of the main goals of a future international linear collider (ILC) is the precise mea- 
surement of the top quark threshold. Next to a precise extraction of the strong coupling 
an unrivaled determination of the top quark mass and its width is possible. This would 
open up a new chapter in the electroweak precision physics which leads to very strong 
checks of the Standard Model or possible extensions. 

The theoretical calculations are based on an effective theory [1,2] (for a review see [3]) 
which enables the simultaneous expansion in the two small parameters present in threshold 
phenomena, the strong coupling a s and the velocity of the heavy quarks, v. Furthermore, 
a resummation of correction terms like (a s /v) n to all orders is possible in perturbation 
theory. In a first step the effective theory is constructed in a matching procedure where 
it is required that the Green functions agree with the corresponding ones in full QCD in 
the limit of infinitely heavy top quark mass. In a second step the Green function in the 
effective theory has to be evaluated leading to the desired cross section. 

Experimental studies [4] have shown that an uncertainty of less than 3% is necessary 
for the theoretical predictions for the cross section of the reaction e + e~ — > ti in the 
threshold region, i.e. for s/s ~ 350 GeV. This has not yet been achieved. 

As far as the strong coupling is concerned currently the analysis is complete to next- 
to-next-to-leading order (NNLO) [5] demonstrating the urgent need for the third-order 
corrections. This huge enterprise has been started already some years ago and many 
building blocks have been evaluated starting from the determination of the third-order 
Hamiltonian [6], third-order corrections to the energy levels and wave functions [7-9] up 
to the evaluation of the contributions from the Coulomb potential to the non-relativistic 
Green function [9]. Next to fixed-order corrections also the resummation of next-to- 
next-to leading logarithmic terms has been considered and quite some progress has been 
achieved [10-13]. 

In order to reach a theory-uncertainty below 3% it is also necessary to worry about 
electroweak effects which can be parametrically as big as the third-order QCD terms. 
First results in this context are available [14, 15] which show that indeed numerical effects 
in the region of a few percent are possible. 

In this paper a further building block needed for the N 3 LO analysis is provided: the 
fermionic corrections to the matching coefficient of the vector current. More precisely we 
compute the three-loop contributions to the matching coefficient containing at least one 
closed loop of massless fermions which we enumerate by n\. 

The paper is organized as follows: in the next section we describe in detail the theo- 
retical framework and the methods used for the calculation. In particular, we classify the 
occurring integrals and discuss their evaluation. Sections El and |U contain a detailed dis- 
cussion about the wave function and matching coefficient up to order a^ni, respectively. 
In Section 03 we briefly discuss the phenomenological impact of our result and Section |H1 
contains our conclusions and outlook. In the Appendix a complete list of master integrals 
is provided — both for the "natural" and the e-finite basis. 



2 



2 Method 



The matching coefficient establishing the relation between the quark current in the full 
and effective theory constitutes a building block for all threshold phenomena involving 
the coupling of a photon or Z boson to heavy quarks. We define the currents through 

f v = QfQ, 

f = 0Vx, (l) 

and generically denote the heavy quark by Q. <fi and x are two-component Pauli spinors for 
quark and anti-quark, respectively. The definition of the matching coefficient is established 
through 

3l = cM? + ^M<P ] ° k D\ + ..., (2) 
6m 2 Q 

where k = 1,2,3 denotes the spacial components. D contains the space-like components of 
the gauge-covariant derivative involving the gluon fields and the ellipsis stand for operators 
of higher mass dimension. The second term on the r.h.s. of Eq. (J2J) is already of NNLO 
with d v = 1 + 0(a s ). The Wilson coefficients c v and d v may be expressed as series in 
a s and represent the contributions from the hard modes which have been integrated out 
in order to arrive at the effective theory. Note that both the higher order terms in the 
inverse heavy quark mass and the time-like component is not of interest for the aim of 
this paper. The purpose of this paper is the evaluation of the fermionic correction to the 
non-singlet contribution of c v containing at least one closed massless quark loop. The 
one- loop corrections to d v contain no contribution proportional to n\. 

In order to compute c v it is obvious to consider the QQ'y vertex which in the following 
is denoted by IV This effectively transforms Eq. (J2J) into 

Z 2 T V = c v Z 2 Z~ l f v + ... , (3) 

where we have on the left- and right-hand side quantities of the full and effective theory, 
respectively, and the ellipses denote terms suppressed by inverse powers of the heavy quark 
mass. In Eq. (jSJ it is understood that the couplings and masses in T v are renormalized. 

Since c v takes into account the degrees of freedom which have been integrated out 
from the full theory it depends in our case only on a s and — for dimensional reasons - 
on In fi 2 /mQ. In particular, c v does not explicitly depend on the momenta of the external 
particles. Thus, it is useful to apply the so-called threshold expansion [16, 17] to Eq. © 
which has the consequence that r„ has to be evaluated for s = ^nig since all except the 
hard region cancel in Eq. Furthermore, on the right-hand side only tree contributions 
have to be considered. In particular we have Z 2 = 1. 

Starting from two-loop order the matching procedure exhibits infra-red divergences 
which are compensated by ultra-violet divergences of the effective theory rendering phys- 
ical quantities finite. In Eq. (J3J) the renormalization constant Z v which generates the 
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Figure 1: Feynman diagrams contributing to the matching coefficient. Bold lines denote 
heavy quarks with mass nig, thin lines denote massless quarks and curly lines denote 
gluons. In (d) and (h) mass counterterm diagrams are shown. 

anomalous dimension of j v takes over this part. Note that the vector current in the full 
theory does not get renormalized. 

Sample Feynman diagrams contributing to the one-, two- and three- loop fermionic part 
of c v are shown in Fig. ^ Due to the special kinematic situation with s = q 2 = and 
on-shell heavy quark lines it is possible to perform a partial fractioning in the integrands 
corresponding to the various types of diagrams and map all occurring integrals to one of 
the following functions: 
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in Ref. [16]. The three-loop functions J±^ and L^' n ^ contain irreducible scalar products 



(2) (2) 

For convenience we have also listed the two- loop functions J± and L± originally defined 
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which are shown as numerators in the three-loop integrals in Eq. (jlj). The corresponding 
indices can only adopt values less or equal to zero. Furthermore, only two out of the three 
indices Tie, n-j and n$ in l}^ n ^ can have positive values. Note that the integrals J+\ lJ+ 
and are actually two-point functions whereas the integrals and Z^' ra '- ) 

correspond to vertices. 

All Feynman diagrams are generated with QGRAF [18]. The various topologies are 
identified with the help of q2e and exp [19,20] which also adapts the notation in order to 
match the one of Eq. (Jlj. 

In a next step the reduction of the various functions to so-called master integrals has 
to be achieved. For this step we have chosen two approaches: the Laporta method [21,22] 
and Baikov's method [23,24] in the formulation of [25]. The application of these methods 
reduces the three-loop integrals to twelve master integrals. Since some of the master 
integrals are only known numerically it is very useful to construct a so-called e-finite 
basis [26]. Details of this procedure and the results for the master integrals are given in 
Appendix El Let us in the following provide more details on our implementation of each 
of the two reduction methods. 

2.1 Implementation of the Laporta algorithm 

The Laporta algorithm is based on the integration-by-parts (IBP) relations [27] where in 
a first step a huge system of equations is created by inserting numerical integer values 
into the relations. After assigning a weight specifying the complexity to each occurring 
integral the system of equations is solved step-by-step. In this way an arbitrary integral 
is expressed in terms of master integrals which cannot be further reduced. 

To date there is only one publicly available computer code, AIR [28], where the Laporta 
algorithm is implemented. We have successfully applied AIR to the two-loop diagrams. 
However, the three-loop integrals of Eq. (J3J) cannot be treated with the help of AIR. 

For this reason we have chosen to write a new implementation of Laporta's algorithm, 
Crusher [29] . It is written in C++ and uses GiNaC [30] for simple manipulations like taking 
derivatives of polynomial quantities. In the practical implementation of the Laporta algo- 
rithm one of the most time-consuming operations is the simplification of the coefficients 
appearing in front of the individual integrals. This task is performed with the help of 
Fermat [31] where a special interface has been used (see Ref. [32]). The main features of 
the implementation are the automated generation of the IBP identities and a complete 
symmetrization of the diagrams. Using Crusher we solved approximately forty million 
equations to reduce all integrals to master integrals. 

2.2 Baikov's method 

Baikov's method is also based on the IBP relations. Here, however, one explicitly con- 
structs the coefficient functions of the master integrals as parametric integrals over some 
polynomial, which encodes the topology of the initial integral. The integrations can then 
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be performed as Cauchy integrations around the origin or as contour integrations be- 
tween the roots of the polynomial. The results are rational functions depending on the 
dimension and the kinematical invariants. 

In some cases, however, it is not possible to perform the integrations in this way. In 
these cases one has to solve the integration-by-parts relations for the parametric integrals. 
They are similar to the relations for the initial integral, however, they depend on less 
indices and are thus in general significantly simpler. Since the relations can still be quite 
complicated, it is convenient to again use Laporta's algorithm. In the current calculation 
the program AIR was used for this task. 

In order to perform the calculations, the algorithms for each integral were implemented 
in FORM [33]. In a first step, however, everything was implemented in Mathematica. This 
is very convenient since it is possible to produce almost all needed Mathematica code 
automatically. Once the program is finished, it can be easily translated to FORM and the 
Mathematica version can be used for debugging. 

The two-loop calculation was done entirely with this method. In the three-loop case, 
however, it turned out that a complete solution of the recurrence relations for the coef- 
ficient functions is not possible with AIR. The calculation was therefore done with the 
Laporta method described above. We used Baikov's method only for the nf part and to 
cross-check some of the coefficient functions. 

3 Wave function renormalization constant 

According to Eq. (jSJ) the wave function renormalization constant in the on-shell scheme 
Z 2 constitutes a crucial input for the computation of c v . It has been computed to two- 
and three- loop approximation in Refs. [34] and [35], respectively. We have repeated the 
calculation of the two-loop and fermionic three-loop contributions which are needed for 
the present calculation and find complete agreement with the literature. For completeness 
we repeat the result for Z 2 which can be cast into the form 
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where L M = ln(/i 2 /m 2 3 ), C F = (iV 2 - 1)/{2N C ), C A = N c , T = 1/2 and n x = n f - 1 is 
the number of light-quark flavours. N c = 3 and a s is the strong coupling renormalized 
in the MS scheme with rif active flavours. £(3) is Riemann's ( function with the value 
C(3) = 1.202057. . ., 7s is Euler's constant with the value 'Je = 0.577216. . . and = 



Li 4 (l/2) = 0.517479 . . .. Note that 5Z K 2 L} and 5Z ( 2 Z> are given to order e 2 and e, respectively, 
which is needed for the three- loop calculation of c v . 

The integrals needed for Z 2 constitute a subset of the ones needed for c v . In fact, from 
Eq. (j3J) only the integral types with a "+" sign contribute. We also want to mention that 
our calculation poses the first independent check of the n\ part of Ref. [35]. 



4 Matching coefficient 

It is convenient to cast the perturbative expansion of the matching coefficient in the form 
c v = l + ^cW+f^Vcf+f-^Vcf + 0(at), (9) 

7T \ 71 / \ 71 J 



8 



where the one- [36] and two-loop [37,38] (see also Ref. [39]) terms are given by 
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with /3q = (llC^/3 — 4Tn//3)/4. Note that the terms proportional to /3oln-*V are con- 

nected to the choice a s (/x) in Eq. (j^J) whereas the ones proportional to tc 2 In originate 

from the separation of hard and soft scales in the construction of NRQCD. 
The decomposition of according to the colour structures is given by 

<4 3) = C F Tni (C pCppL + C A c FAL + Tc FHL + TniC FLL ) + non-n, terms . (11) 

In the notation of Eq. © it is understood that in r„ the renormalization of the cou- 
pling constant and heavy quark mass is already performed. The renormalization of the 
coupling constant can be done straightforwardly by replacing the bare coupling with the 
renormalization constant times the renormalized coupling. The mass renormalization, 
however, is more intricate since we are dealing with on-shell integrals. It is therefore con- 
venient to calculate the counterterms directly by considering one- and two-loop diagrams 
with zero-momentum insertions in the massive fermion lines. Sample diagrams are shown 
in Fig. dd) and (h). The vertex denoted by a cross is replaced by the mass renormal- 
ization constant in the on-shell scheme. One-loop diagrams with two insertions are not 
needed at the order considered in this paper. 

Starting from three-loop level the wave function renormalization constant begins to 
depend on the QCD gauge parameter, £ [35]. However, it is interesting to note that the 
colour structures entering our result (cf. Eq. (jHJ)) are independent of £. As a consequence 
the quantity T v also has to be independent of £ which is indeed the case in our calcula- 
tion. Actually, the £ dependence remaining in the sum of all genuine three-loop diagrams 
contributing to r„ is canceled by the contribution where a two-loop mass counterterm is 
inserted in the one- loop vertex diagram (cf. Fig. IHd)). Note that the mass counterterm 
contribution from the two-loop fermionic diagram (cf. Fig. ^h)) is £ independent. The 
cancellation of the gauge parameter serves as a welcome check for the correctness of our 
result. 

As already mentioned above, after taking into account all counterterm contributions 
there are still infra-red divergences left in the quantity Z 2 T V and thus in c v Z~ l . This is a 
consequence of the threshold expansion accompanied by dimensional regularization which 
is used to extract the matching coefficient. Alternatively it would have been possible to 
choose a cut-off for the momentum integrations which finally results in a factorization scale 



9 



separating the hard and soft momenta. In our approach the infra-red poles are absorbed 
into Z v which we define in the MS scheme. As a consequence the matching coefficient 
itself is finite but scale dependent. In the physical quantities this scale dependence gets 
canceled against the corresponding contributions from the effective theory. Of course, one 
has to make sure that for the subtraction of the poles in the effective theory the same 
scheme is used as in the full theory. 

For the renormalization constant Z v we obtain 
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where the ellipses stand for non-n^ and O {a A s ) terms. 
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The uncertainties assigned to the numerical constants in c FFF and c F al are based on a 
conservative estimate. Note that the precision of these quantities is more than enough for 
all phenomenological applications. Inserting the numerical values for the colour factors 
we obtain for \i = mq 
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The coefficients in Eq. (|13|) correspond to an expansion parameter a s (/i), as given in 
Eq. (0). Choosing instead « s (mg) leads to 
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The dependence on L M in Eq. (JT5|l is canceled against contributions from the effective 
theory. They agree with the ones of Ref. [40] 1 . 



5 Phenomenological applications 



Let us in this section estimate the numerical effect of our new terms on the bottom and 
top system. In particular we consider the decay of the T(IS') bound state to leptons and 
the production of top quark pairs close to threshold. 

Next to the matching coefficient considered in the previous Sections a crucial ingredient 
for these quantities is the wave function at the origin. Currently the second order is known 
completely [41-45] and at order the quadratically [46,47] and linearly [12,40] enhanced 
logarithms and the corrections proportional to /3q [8,9] are available. It is convenient to 
introduce the quantity p\ = \ipi(0)\ 2 /\ipf (0)\ 2 where the Coulomb wave function is given 
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where fi s = C F m q a s (n s ) is the soft scale. It is straightforward to obtain the result for 
general a s (/x) using standard renormalization group analyses. In Eqs. (fIH jl -(fI5 )l S is the 
spin quantum number which is equal to one in our applications. The ellipses in Eq. (II 6 j) 
represent yet unknown corrections like, e.g., the pure ultrasoft contributions. 

For completeness we provide the one- and two-loop coefficients of the (3 function and 
the static potential which are given by 
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5.1 Bottom system 

The leptonic decay of the T(IS') state can be cast in the form [41-45] 
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with T LO = 4vriV c Qga 2 |<0f(O)|7 (3ml), Qb = -1/3, N c = 3, a is Sommerfeld's fine- 
structure constant and nonperturbative contributions to Eq. (}2*T]) are ignored. 
Inserting the perturbative expansion for p\ and c v we obtain 

Ti » T\° (l - 1.70 a s {m b ) - 7.98 a 2 s {m b ) + 30.0 a 3 s {m b )\ ni + .. .) 
x 1 - 0.30 a a (Vs) + ol 2 s {hs) (17.2 - 5.19 In a,(ji a )) 
WM (-14-4 In 2 a a (jji s ) + 0.17 In a s (/x s ) - 34.9^) + . . .] , (22) 

where in the matching coefficient \x = m b has been chosen and the corresponding strong 
coupling is defined with five and a s (/i s ) with four active flavours. 

Starting from a s (Mz) = 0.118 we have used the program RunDec [48] to obtain 
a s (m b ) = 0.2096 and a,(ji a ) = 0.2967 with m b = 5.3 GeV and fi s = 2.0967 GeV which 
leads to 

ri « r^°(l - 0.446nlo + 1.75nnlo - 1.20 N3LO / + ...), (23) 

where Eq. f!22j) is expanded and terms of order are dropped consistently. The prime re- 
minds that the third-order corrections are not complete. Apart form the new contribution 
to c v and the known third-order corrections to p\ we have also included all interference 
terms which are proportional to powers of n\. Note that the new corrections are responsi- 
ble for the reduction of N 3 LO' terms from —1.67 to —1.20 which amounts to about 47% of 
the Born cross section. However, in total the fermionic corrections tend to further reduce 
the strong increase of the perturbative coefficients leading to an overall correction factor 
of approximately 10%. 

We have performed the numerical analysis also for m b = 4.8 GeV. There are changes 
in the individual contributions to I\ (c.f. Eq. (J23|) ) of the order of a few per cent, however, 
the final correction factor remains the same. 



5.2 Top system 

In the top-quark case, the nonperturbative effects are negligible. However, the effect of 
the top-quark total decay width T t has to be properly taken into account [49], as it is 
relatively large and smears out the Coulomb-like resonances below threshold. The NNLO 
analysis of the cross section [5] shows that only the ground-state pole gives rise to a 
prominent resonance. 

A crucial quantity in connection to the threshold production of top quark pairs is 
the peak of the normalized cross section R = <r(e + e~ — > tt)/cr(e + e~ — > It is 

dominated by the contribution from the would-be toponium ground-state, which can be 
cast in the from 



i?i(e + e" 



tt) 



K l Pi 



c 2 v {m t ) + 



C 2 F a 2 s (n t 
12 



-c v (m t ) (d v (m t ) + 3) 



+ . . . , (24) 
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with the leading order term R\° = 6iiN c Q^\ , ipf (0)\ 2 / {m 2 Y t ). The contributions from the 
higher Coulomb-like poles and the continuum are not included in Eq. 
The analog equations to Eqs. (|22|) and (}23j) read 

R 1 » (1 - 1.70 a s (m t ) - 7.89 a 2 s (m t ) +37.2 a 3 a (m t )\ ni 

1 - 0.43a s (/i s ) + o?M (16.1 - 5.191na.(^)) 



+ ...) 



x 



+a 3 M (-13.81n 2 a s (/i s ) +2.06 \na s (fi s ) -27.2|^J + 



(25) 



where a s (m t ) and a s (/j, s ) are defined with six and five active flavours, respectively. Using 
again a s (Mz) = 0.118 one gets a s {m t ) = 0.1075 and a s (/j, s ) = 0.1398 with m t = 175 GeV 
and /j, s = 32.625 GeV and finally 



Ri 



- 0.243 NLO + 0.435 N nlo - 0.195 N 3 Lcy + 



(26) 



The fermionic corrections to c v are responsible for a reduction of the third-order coefficient 
from —0.268 to —0.195 and thus amount to moderate 7% of the leading order term. 
Similarly as for the bottom quark case also for the top quark the perturbative series is 
alternating and the third-order coefficient tends to stabilize the expansion. It is interesting 
to note that after the inclusion of our new terms the total corrections amount to less than 
1%. 



6 Conclusions and outlook 

This paper deals with the question of establishing a relation for the vector current between 
full QCD and NRQCD. The corresponding matching coefficient, c v , constitutes a building 
block in all threshold phenomena involving a vector coupling. The main result of this 
paper is the fermionic three-loop contribution to c v which contains a closed light fermion 
loop. The numerical effect of the new terms is relatively big and amounts to about 47% 
for the bottom and to about 7% for the top quark system. However, their inclusion leads 
to a reduction of the overall corrections. E.g., in the case of the top quark the corrections 
to the normalization at the peak of the production cross section amount to less than 1% 
after including all currently known perturbative terms. 

Several steps still have to be taken in order to arrive at the complete analysis of the 
top-anti-top quark system at threshold. Next to the non-fermionic pieces to c v the major 
building blocks which are still missing are the ultra-soft corrections and the three-loop 
static potential. 
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Figure 2: Three- loop master integrals. Bold lines denote massive lines with mass 2mg, 
thin lines denote massive lines with mass rriQ and dashed lines denote massless lines. All 
external lines are on-shell. A dot on a line denotes a squared propagator. 

A Results for the master integrals 

The master integrals for the three-loop functions of Eq. are shown in Fig. El For the 
self-energy integrals the external momentum squared equals the square of the heavy 
quark mass, q 2 = nig. For the integrals I7 and /§ we have q 2 = 4:mQ. The momenta 
qi and q^ flowing out of the vertex diagrams on the right fulfill q\ = q\ 
(<7i + Q2) 2 — ^ m Q- 111 Minkowski space the results for the master integrals read 

h = jf(0, 1,0, 1,1, 0,1, 0,1) 



tuq with 



3c 



in 



+ 



302 



2tt 2 



1 

3e 3 
89 



3 
4' 



j. / 3059 , 

— C(3) H 7T 4 1 e 

3 SV ; 1440 



-734 



10 11 
~3~ ' ! 

69 



25 

y 



-7T 



16C(3) 



9751 
1440' 



-7T 



107 2309 wr . 
+— tt 2 C(3) + -— C(5) 



4 



+ O (e 3 ) 



jf(0, 1,1, 1,0, 1,1, 0,0) 



3c 



^ 2 \ T 4 (l - e)r 2 (e) r(3e - l)r(3 - 6e) 



m 



Q , 



T 2 (2 - 2e) 



T(3 - 4e) 



(27) 



(2f 



15 



Lf nj) (0,0, 1,1, 1,1, 0,1,0) 
2 ( V 2 ^ \ 3e f 2 10 /26 1 ,\ 1 9 2 14 , , 

(£| e J {3? + 3^ + (t + T ) ! + 2+ 2* + T c(3) 

/ 398 53 2 2, 238 w x 287 A / 259 2 

+ ( 3~ + ~2~ ~ ^ M 2 + " C(3) + 720^ J £ + (" 1038 + IT 

-160. 2 In 2 + ^ In 2 2 + ^ C (3) + fj^ + ^ 2 C(3) + ^C(5) 



+ yln 4 2 + 512a 4 ^) e 2 + (e 3 ) J , 
jf (0,0,0,1,1,0,1,0,1) 



Jf (0,0,1,1,1,0,1,0,0) 



m « I &f 1E ) 3 { h + h + (f + H 7 + 33 - 8328(4) + 152 - 870(4)e 



(29) 



4 fV \ 3e r 2 (i- e )r(6)r 2 ( 2e -i)r( 3e -2) 

r(4e-2)r(2- £ ) ' (30) 

j| 3) (0,l,0,l,0,0,l,0,l) 

m 4 ( f X\ (c r 2 (i- g )r( g )r( 2g -i)r(3- 4e ) 

Q V"?/ r(2-2e) r(3 - 3e) ' (31) 
r(3) 



/i 2 \ r 3 (i- e )r(26-i)r(36-2)r(5-66) 

Q V^q7 r(3-3e)r(4-4e) ' ^ 

Ji 3) (0,0,0,l,l,0,l,0,l) 

m 4 [ - 7B ] f 1 + 1 + /_1Z 1 2 \ 1 _ 6 . 68 273 8 7(i) 

V ™g / 1 3e 3 2e 2 V 36 12 / e 

-56.300353(l)e - 209.48231(l)e 2 + O (e 3 ) J , (33) 
J^O, 0,0, 2, 1,0, 1,0,1) 

-^(^•) 3, {i+ 3 l + (-I + H^ ia797602(i) 

+62.250613(l)e + (e 2 ) J , (34) 
Ji 3) (0,l,0,l,l,0,l,0,l) 



+0 (e 2 ) } , (35) 
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£ 10 



4 3 ' n,) (o,o,i,i,i,o,i,i,o) 



= < (£r 7b ) + & + ( 8 + H -e + 52 - 5698(4) + 145 - 087(4)e 

+562.250(14)e 2 + O (e 3 ) | , (36) 

r n = ji 3) (o,i,o,i,o,i,i,o,i) 

1 29 \ 1 

— + — 7T 2 - + 34.791(4) + 154.08(2)e 
6 24 / e 

+0 (e 2 ) } , (37) 

ri2 = / 3) (o, 1,0, 1,0, 1,1, 0,2) 




£re-J + ^ + Q + ; + 7.024(4) + 32.04(2), 

+0 (e 2 ) } , (38) 

where C(5) = 1.036927 . . . and a 4 = Li 4 (l/2) = 0.517479 . . .. 

Ii was calculated in Ref. [35]. We have repeated the calculation described in the 
reference with the program XSummer [51] and find complete agreement. The results for 
/3-/6 can be found in Ref. [21]. We have checked the result for ^3 numerically with 
the help of a one-fold Mellin-Barnes [52,53] representation and find complete agreement. 
Note, that l2,h, h an d Iq are quite simple to evaluate and are available for general e. 
Nevertheless the result for J 6 as given in Ref. [35] is not correct whereas the correct result 
can be found in Ref. [21]. The pole part of J7 agrees with the result of Ref. [54]. 

The remaining master integrals, I7-I12, have been evaluated with the help of the 
Mellin-Barnes method where the evaluation of the integrals has been performed with 
the program MB [55]. I 7 and I 8 can be expressed in terms of a one-fold Mellin-Barnes 
representation and are thus known with a quite high precision. On the other hand, iio is 
represented by a two-fold and J9, In and I12 even by a three- fold integration which results 
in less accurate results. The quoted uncertainties in the above equations correspond to 
twice the Vegas error given by MB for the multi-dimensional integrals and to a conservative 
estimate in case of the one- dimensional integrals I7 and Is- Note that the latter errors 
are negligible. 

As can be seen in the above results some of the master integrals are needed to higher 
order in e which on one hand makes the calculation very tedious and on the other hand 
leads to less accurate results. For this reason we decided to change basis and switch - 
at least for some of the integrals — to the so-called e-finite master integrals [26] which 
have the advantage that the coefficient function is finite and thus the integral itself is only 
needed to order e°. 
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Figure 3: e-finite master integrals. The same coding as in Fig. [21 is adopted. 



Since I\-I& are known analytically, a replacement is, of course, not necessary. Fur- 
thermore, for I 7 and Ig the numerical precision is sufficient for our calculation. As far 
as the remaining four integrals are concerned, we found it convenient to replace Ig and 
1 12 by the integrals shown in Fig. El Their numerical evaluation with the help of MB is 
straightforward leading to the results 

4 = J® (0,1, 1,1, 1,0, 1,0,1) 

- (^) 3t {^ + i + (f + H7 +64 - 678(8)+o(f) }- (39) 

l[ 2 = J W (0, 1,0, 1,1, 1,1, 0,1)= (Ke-A | 2 ' 4 °f W + 8.1(2) + 0(e)j .(40) 

f f 

Switching from Ig and I12 to Ig and I l2 reduces the number of coefficients which are 
only known numerically from 17 to 14. In particular the 1/e-poles of Ig and I\2 have 
been determined in analytical form from the e-finite integrals via the corresponding IBP 
relations. 

A further reduction is achieved after exploiting the locality and incorporating the 
knowledge of the linear ln(/i 2 /mg) term in c v [40] which we have checked numerically in 
a first step. In this way the analytical result for the 1/e-pole of In has been determined. 
In total ten (eleven) numerical coefficients contribute to our final result given in Eqs. (|13|) 
and (fTlj) using the e-finite ( "normal" ) basis. 

In the end it turned out, that the numerical precision is better using the "normal" 
basis. The reason is, that the uncertainty of the e°-coefficient of l{ 2 is by far the largest. 
Still, it was useful to consider the e-finite basis since it enabled us to determine analytical 
results for two integral coefficients of the "normal" basis. 
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